Mean Sea Surface Temperature#

https://www.ncei.noaa.gov/products/optimum-interpolation-sst

Hide code cell source
import warnings
warnings.filterwarnings("ignore")

import sys
import os
import os.path as op

import xarray as xr
import plotly.express as px
import geopandas as gpd
import pandas as pd
import cartopy.crs as ccrs
import numpy as np
import matplotlib.pyplot as plt

sys.path.append("../../../../indicators_setup")
from ind_setup.plotting import plot_base_map, plot_bar_probs, plot_map_subplots, add_oni_cat
from ind_setup.plotting_int import plot_timeseries_interactive, plot_oni_index_th

sys.path.append("../../../functions")
from data_downloaders import download_oni_index
lon_site, lat_site = 134.368203, 7.322074

#Area of interest
lon_range  = [129.4088, 137.0541]
lat_range = [1.5214, 11.6587]
shp_f = op.join(os.getcwd(), '..', '..','..', 'data/Palau_EEZ/pw_eez_pol_april2022.shp')
shp_eez = gpd.read_file(shp_f)
path_data = "../../../data"
path_figs = "../../../matrix_cc/figures"
data = xr.load_dataset(op.join(path_data, 'sst_daily_1981_2024_palau.nc'))
dataset_id = 'sst'
ax = plot_base_map(shp_eez = shp_eez, figsize = [10, 6])
im = ax.pcolor(data.lon, data.lat, data.mean(dim='time')[dataset_id], transform=ccrs.PlateCarree(), 
                cmap = 'hot_r', vmin = np.percentile(data.mean(dim = 'time')[dataset_id], 1), 
                vmax = np.percentile(data.mean(dim = 'time')[dataset_id], 99))
ax.set_extent([lon_range[0], lon_range[1], lat_range[0], lat_range[1]], crs=ccrs.PlateCarree())
plt.colorbar(im, ax=ax, label='SST (ºC)')

plt.savefig(op.join(path_figs, 'F12_SST_map.png'), dpi=300, bbox_inches='tight')
../../../_images/d5e2a63926bbfe0c9d8d451cea54d8ea32aec49f38d4a415f62a1995c8580e83.png

Annual average#

data_y = data.resample(time='1YE').mean()
im = plot_map_subplots(data_y, dataset_id, shp_eez = shp_eez, cmap = 'hot_r', 
                  vmin = np.nanpercentile(data_y.min(dim = 'time')[dataset_id], 1), 
                  vmax = np.nanpercentile(data_y.max(dim = 'time')[dataset_id], 99),
                  cbar = 1, cbar_pad = 0.05)
../../../_images/8ecca7714035b3c56d8a3efa9234cad327a296a0ca03e84db8c44732652f0011.png

Annual anomaly#

data_an = data_y - data.mean(dim='time')
fig = plot_map_subplots(data_an, dataset_id, shp_eez = shp_eez, cmap='RdBu_r',
                        vmin=-2, vmax=2, cbar = 1, cbar_pad = 0.05)
../../../_images/865ee4492f852129e25cbf141ca2dcd03cb799e2d92494d9390535563bc0a6c4.png

Average Area#

data_mean = data.mean(dim = ['lon', 'lat']).to_dataframe()
datag = data_mean.groupby(data_mean.index.year).max()
datag.index = pd.to_datetime(datag.index, format = '%Y')
datag['sst_ref'] = datag['sst'] - datag.loc['1961':'1990'].sst.mean()
data_mean = data.mean(dim = ['lon', 'lat']).to_dataframe()
data_mean = data_mean.resample('YE').mean()
data_mean['sst_ref'] = data_mean['sst'] - data_mean.loc['1961':'1990'].sst.mean()
fig, ax, trend = plot_bar_probs(x=data_mean.index.year, y=data_mean.sst_ref, trendline=True,
                                y_label='SST [°C]', figsize=[15, 4], return_trend=True)
ax.set_title('Annual Mean SST anomalies (Over and above 1961 - 1990 reference period)', fontsize=15);
plt.savefig(op.join(path_figs, 'F12_SST_trends_Annomalies.png'), dpi=300, bbox_inches='tight')
../../../_images/a5325ab3fbee3ceb082ca461e0c5db7753853f6255f6ba6be35ecc72ef502284.png

Given point#

loc = [7.37, 134.7]
dict_plot = [{'data' : data.sel(lon=loc[1], lat=loc[0], method='nearest').to_dataframe(), 
              'var' : dataset_id, 'ax' : 1, 'label' : 'SST [ºC]'},]
ax = plot_base_map(shp_eez = shp_eez, figsize = [10, 6])
ax.set_extent([lon_range[0], lon_range[1], lat_range[0], lat_range[1]], crs=ccrs.PlateCarree())
ax.plot(loc[1], loc[0], '*', markersize = 12, color = 'royalblue', transform=ccrs.PlateCarree(), label = 'Location Analysis')
ax.legend()
<matplotlib.legend.Legend at 0x183943110>
../../../_images/038a2757803511e8b284fa11a57d0ba084269dbabe56e09215299e69572e231d.png
fig = plot_timeseries_interactive(dict_plot, trendline=True, scatter_dict = None, figsize = (25, 12), label_yaxes = 'SST [ºC]');

ONI index analysis#

p_data = 'https://psl.noaa.gov/data/correlation/oni.data'
df1 = download_oni_index(p_data)
lims = [-.5, .5]
plot_oni_index_th(df1, lims = lims)

Group by ONI category

data_xr = xr.load_dataset(op.join(path_data, 'sst_daily_1981_2024_palau.nc')).rename({'lat':'latitude', 'lon':'longitude'}).resample(time='1M').mean()
data_xr['time'] = data_xr.time.values.astype('datetime64[M]')

data_xr['ONI'] = (('time'), df1.iloc[np.intersect1d(data_xr.time, df1.index, return_indices=True)[2]].ONI.values)
data_xr['ONI_cat'] = (('time'), np.where(data_xr.ONI < lims[0], -1, np.where(data_xr.ONI > lims[1], 1, 0)))
data_oni = data_xr.groupby('ONI_cat').mean()
data_oni['labels'] = (('ONI_cat'),['La Niña', 'Neutral', 'El Niño'])
dataset_id = 'sst'
fig = plot_map_subplots(data_oni, dataset_id, shp_eez = shp_eez, cmap = 'hot_r', 
                  vmin = np.percentile(data_xr.mean(dim = 'time')[dataset_id], 0) - .25, 
                  vmax = np.percentile(data_xr.mean(dim = 'time')[dataset_id], 100) + .25,
                  sub_plot= [1, 3], figsize = (15, 8), cbar = True, cbar_pad = .1)

plt.savefig(op.join(path_figs, 'F14_SST_ENSO.png'), dpi=300, bbox_inches='tight')
../../../_images/e3f56ecbd4b2cada7c4a664bd4c1b6e03ee614733842c3ba96ad845abf3c67e2.png
data_an = data_oni - data_xr.mean(dim='time')
data_an['labels'] = (('ONI_cat'),['La Niña', 'Neutral', 'El Niño'])
fig = plot_map_subplots(data_an, dataset_id, shp_eez = shp_eez, cmap='RdBu_r', vmin=-.4, vmax=.4,
                  sub_plot= [1, 3], figsize = (15, 8), cbar = True, cbar_pad = 0.1)
../../../_images/6ae2c6bcb79a79ce48d7b6e195bef095ae22793262242eabed64c85077b97424.png